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ABSTRACT 

Context. Coronal loops act as a resonant cavity for low frequency fluctuations that are transmitted from the deeper layers of the solar 
atmosphere. Such fluctuations are amplified in the corona and lead to the development of turbulence that in turn is able to dissipate 
the accumulated energy, thus heating the corona. However trapping is not perfect, some energy leaks down to the chromosphere on a 
long timescale, thus limiting the turbulent heating. 

Aims. We consider the combined effects of turbulence and energy leakage from the corona to the photosphere in determining the 
turbulent energy level and associated heating rate in models of coronal loops which include the chromosphere and transition region. 
Methods. We use a piece-wise constant model for the Alfven speed in loops and a Reduced MHD - Shell model to describe the 
interplay between turbulent dynamics in the direction perpendicular to the mean field and propagation along the field. Turbulence 
is sustained by incoming fluctuations which are equivalent, in the line-tied case, to forcing by the photospheric shear flows. While 
varying the turbulence strength, we compare systematically the average coronal energy level and dissipation in three models with 
increasing complexity: the classical closed model, the open corona, and the open corona including chromosphere (or 3-layer model), 
the latter two models allowing energy leakage. 

Results. We find that: (i) Leakage always plays a role: even at for strong turbulence, the dissipation time never becomes much lower 
than the leakage time, at least in the three-layer model. Hence, the energy as well as the dissipation levels are systematically lower 
than in the line-tied model, (ii) In all models, the energy level is close to the resonant prediction, i.e., assuming effective turbulent 
correlation time longer than the Alfven coronal crossing time, (iii) The heating rate is close to the value given the ratio of photospheric 
energy divided by the Alfven crossing time, (iv) The coronal spectral range is divided in two, an inertial range with 5/3 spectral slope, 
and a large scale peak where nonlinear couplings are inhibited by trapped resonant modes, (v) In the realistic 3-layer model, the 
two-component spectrum leads to a global decrease of damping equal to Kolmogorov damping reduced by a factor u„„JV^ where 
is coronal Alfven speed. 

Key words. Sun: corona, transition region - Magnetohydrodynamics (MHD) - Turbulence - waves - Methods: numerical 

' 1. Introduction incompressible, so the transverse structure in density may be ne- 

glected compared to the gravitational stratification, while the ex- 

^ Solving the coronal heating problem involves understanding ^^^^^^^ ^^e field from the denser layers of the photosphere 

d how fast magnetic energy can be accumulated in the corona and chromosphere into the corona may be taken into account 

how fast this energy is dissipated. We investigate this problem gradients along the field of the Alfven speed. The resulting, 

by considering a model loop in which kinetic and magnetic en- simplified coronal loop retains the basic ingredients which lead 

ergies are injected into the corona in the form of Alfven waves heating: turbulent coupling and propagation through a strat- 

generated by photospheric motions. A large body of work has ^^^^ atmosphere where stratification appears as an increase of 

been devoted to this problem, (Milano et al. 1997; Dmitruk et al. ^jf^g^ ^pg^^ ^^^^ photosphere to corona. 

2003; Rappazzoetal. 2007, 2008; Nigroetal. 2004, 2005, 2008; ^he stratification is characterized by the ratio of mean 

Buchhn & Velli 2007): we consider here a previously neglected ^jfy^^ ^^^^^^ i„ ^^e photosphere (V^) and in the corona (VD 

effect which plays a large role in regulating the turbulent energy ^j^i^j^ i^ ^ ^j^^^jj parameter: 

balance in the corona, namely the leakage of coronal energy back ^ ^ 

down to the photosphere. ^ - ^a/K « ^ (1) 

A solar loop can be described as a bundle of magnetic field The part of the wave spectrum incoming into the corona that 

lines that expand into the corona but are rooted in the denser we shall consider here is the low frequency part, for which the 

photosphere at two (distant) points, so that their length is typi- Alfven speed contrast is seen by waves of frequency w as a sharp 

cally much greater than the transverse scale. The magnetic field transition. This occurs if 
is therefore mostly along the direction of the loop, and provided 
the transverse magnetic field is not too strong the curvature of the 

loop may be neglected. In addition, if the ratio of the plasma to for ~ 2000 km/ s and a transition region thickness of about 

magnetic field pressures is small the motions are predominantly H = 200 km. For these low frequencies, the transition region 



cj < max(IVVal) * (V', - V")/// * 5 - 10 Hz (2) 
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(T.R.) acts as a transmitting and reflecting barrier, with the im- 
portant property that the transmission is not symmetric, so that 
a coronal loop acts as a cavity which resonates at specific fre- 
quencies, based on the Alfven crossing time fjj - LdV^ (Lc is 
the length of the coronal part of loop): 

OJ = nnV^/Lc = nn/f^i (3) 

with n = 0, 1... (lonson 1982; Hollweg 1984). 

The cavity is perfectly insulated in the limit of infinite Alfven 
speed contrast, i.e. e - 0, which corresponds to the so-called 
line-tied limit. In this limit, the corona exerts no feedback on the 
solar surface. The zero frequency resonance is clearly distinct 
from the finite frequency resonances; in the former, the coro- 
nal magnetic energy grows without bounds while the kinetic en- 
ergy remains finite (Parker 1972; Rappazzo et al. 2007); in the 
latter case, both magnetic and kinetic coronal energies grow at 
equipartition. 

In reality, the trapped energy is limited, because the cavity 
looses energy by two different mechanisms: damping (turbulent 
or not), and leakage, due to the finite Alfven speed contrast. The 
leakage time is given by (Hollweg 1984; Ofman 2002; Grappin 
et al. 2008): 

tL = LjV^ (4) 

Note that the leakage time is much greater than the Alfven cross- 
ing time, since = eti ' . The dissipation rate of the loop will 
thus depend on (i) the energy input into the corona, as well as 
its frequency distribution (resonant or not) (ii) what part of the 
energy input goes into heat what part returns back to the solar 
surface (leakage). 

In the previous works starting with Hollweg (1984), it has 
always been assumed that the leakage time was long compared 
to the (turbulent) dissipation time, thus leakage was neglected 
(line-tied limit). Because neglecting leakage implies neglecting 
the back-reaction of the corona on the deeper layers, in the line- 
tied limit the velocity can be imposed at the coronal base. This 
is justified if the leakage time is larger than the coronal dissipa- 
tion time. Estimating the latter to be given by the photospheric 
turnover time fj^^ = Ij_/{27tUo), we have for the ratio of the two 
time scales: 

XL = fL/f^L - iL/L)2MUo/V^) > 2n (5) 

Since the coronal energy per unit mass is expected to reach larger 
values than at the surface, this largely justifies neglecting leak- 
age. However, identifying the dissipation time with the turnover 
time might be erroneous, as turbulence, at least in some simula- 
tions (e.g., Nigro et al. 2008) shows a high degree of intermit- 
tency, so that the dissipation time is orders of magnitude larger 
than such simple estimates. 

This motivates us to relax the line-tied hypothesis, using 
models of turbulent loops that include leakage. The problem be- 
comes then more complex, as the velocity boundary conditions 
are no longer fixed, the velocity being the sum of the incom- 
ing coronal base field and the outcoming coronal signal. We will 
consider two versions of the problem that includes leakage. In 
the first version, which will be called the one-layer model, we 
simply change the boundary conditions at the coronal base, tak- 
ing leakage into account. The incoming spectrum depends partly 

' As will be seen in sec. 2.3.3, eq. 25, at every reflection a fraction e 
of the coronal energy leakes from the transition region down to the chro- 
mosphere, so one needs 1/e reflections to evacuate the coronal energy, 
i.e. a timescale fjj/e. 




(a) line-tied model (b) 1/2 open corona (c) 1/2 open corona 

(closed) (one-layer) (three-layer) 



Fig. 1. Sketch of the coronal heating process. Above : the general prob- 
lem of photospheric injection, transmission, turbulent dissipation, and 
leakage back to the photosphere. Red thick arrows at the left foot point 
representing the surface shear forcing. Below: the three numerical mod- 
els considered in this paper: (a) closed model (no leakage) with imposed 
velocity at the coronal base (b) semi-transparent corona with imposed 
wave input at the coronal base (c) semi-transparent corona including 
chromospheric turbulence, with imposed wave input at the chromo- 
spheric base. Thin arrows indicate the wave reflection and transmission, 
white thick arrows represent the leakage out of the numerical domain. 



of the (given) signal assumed given by the chromospheric lay- 
ers below, partly on the signal propagating downward from the 
corona and being largely (but not fully) reflected. In the second 
version, which will be called the three-layer model, the domain 
is enlarged to include two chromospheric layers. In that case, the 
signal propagating upward from the coronal base is still more 
uncontrolled than in the previous case, as the chromospheric tur- 
bulence which develops and determines the state of the coro- 
nal base is not directly predictable from the photospheric input. 
Fig. 1 summarizes the models: the classical closed model, and 
the two versions including leakage. 

To describe the turbulence dynamics along the loop, we 
will use the Shell model for Reduced MHD (Nigro et al. 2005; 
Buchlin & Velli 2007). Shell models of turbulence share with 
full turbulence power-law energy spectra, as well as chaotic (in- 
termittency) properties which are very close to direct numerical 
simulations of primitive MHD equations (Gloaguen et al. 1985; 
Biskamp 1994) The system will be forced by introducing DC 
fluctuations, i.e., a spectrum of fluctuations at different perpen- 
dicular scales that is constant in time. 

We will show that the finite leakage time leads to signifi- 
cant differences with previous results obtained using line-tied 
boundary conditions. The plan is the following. The next sec- 
tion deals with basic physics, model equations and parameters. 
Section three deals with simple phenomenology. Results are in 
section four, section five contains the discussion. 

2. Basic physics, model equations and parameters 

2.1. Three-layer atmosphere: linear reflection/transmission 
laws 

We begin by describing our model atmosphere and the prop- 
erties of linear Alfven wave propagation within such an at- 
mosphere. The atmosphere is considered to be stratified in the 
vertical direction, with three successive layers representing a 
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Incident 
signal 
Reflected 
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Fig. 2. The three-layer model: sketch of the transmission and reflection 
properties of transverse fluctuations at the coronal bases of a magnetic 
loop with piece-wise constant Alfven speed, in the particular case con- 
sidered here (no input from right chromosphere). 



left photosphere/chromosphere, the corona, and a right photo- 
sphere/chromosphere. The atmosphere is threaded by a vertical 
uniform field Bo along which Alfven waves propagate. In each 
of these three layers, the Alfven speed is constant, so that a pro- 
gressive Alfven wave propagates at constant speed without de- 
formation. When a wave encounters a density jump interface, the 
velocity and magnetic field fluctuations, which are parallel to the 
interface, are continuous. The proper Alfven modes propagating 
in opposite directions along the loop are defined by the Elsasser 
variables: 



- U + bj yfp 



(6) 



where p is the density and u,b aie the velocity and magnetic 
field fluctuations, which are in planes parallel to the photo- 
sphere/corona transition region. Assuming a positive mean field 
Bq, the quantity will propagate to the right and the quantity 
z~ to the left. It is immediately seen from this definition that the 
density jump at the transition region will determine a wave am- 
plitude jump of order Ij ^ - 1/e. The derivation of the jump 
relations may be found in (Hollweg 1984). Continuity of the ve- 
locity and magnetic field fluctuations at the two interfaces imply 
the following relations between wave amplitudes respectively at 
left and right boundaries: 



•^1 



(4 



ZR)e 



(7) 



We use 1, L to denote the amplitudes at the left T.R. (resp. 1 on 
the photospheric side, L on the coronal side), and 2>,R to denote 
the amplitudes at the right T.R. (resp. R on the coronal side, 3 on 
the photospheric side); see Fig. 2: 



zf = z±(x = 0-) 
zl = z±(x = 0+) 
zt = z\x = L-) 
zt = z\x = L^) 



(8) 
(9) 
(10) 
(11) 



with the exponents 4- or - in and indicating whether we are 
on the right or the left side of the two transition regions, located 
respectively ?X x - Q and x - Lc- 

To obtain the jump conditions to be effectively implemented 
in the three-layer model, we rewrite Eqs. 7 as follows. We de- 
note by input what goes into the corona and output what goes 



out. The coronal inputs z^ and z^ are expressed in terms of the 
chromospheric inputs (z| and zp and the coronal outputs (z^ 
and zp. Similarly the reflected chromospheric signals z\ and Z3 
are expressed in terms of the chromospheric inputs and of the 
coronal outputs: 

z^ = (1 -H a)z\ - azi 

= (1 -H fl)Z3 - flZg 

Zj" = (1 - a)zi + flz| 
Z3 = (1 - a)z^ + aZ3 

The parameter a: 

a = (1 - 6)/(l + e) 



(12) 



(13) 



is the reflection coefficient. It is instructive to consider the limit 
e = 0. Then the coronal reflection coefficient a becomes unity. In 
this case, the velocity at the left coronal boundary is exactly z\, 
that is, specifying the chromospheric input is the same as speci- 
fying the velocity (and the same at the right coronal boundary): 
this is the well-known line-tied limit. In this limit, the magnetic 
field fluctuation is not specified and depends on the coronal evo- 
lution, since one has: bil ^[p - -z|+z^. Returning to the general 
case with a non-zero Alfven speed ratio e, we see that specifying 
the chromospheric input does not directly determine the velocity 
at the T.R. either. We will choose here to consider a non-zero 
input only from the left foot point (boundary), in order to better 
follow the propagation of the incident signal. 

In the early work by Hollweg (1984), the three-layer model 
was studied analytically, with a damping term representing the 
effects of turbulence. As we said, turbulent dissipation is highly 
intermittent thus requiring a description that goes beyond a sim- 
ple damping term. We now define the nonlinear part of the 
model, i.e., the turbulence model. 

The jump conditions just described are not specific of a lin- 
ear framework. In the general case where the waves have a per- 
pendicular structure and interact nonlinearly, the jump condi- 
tions hold as well. In the final model to be detailed now, where 
the wave amplitudes depend on the coordinate along the loop 
and on an index n representing the perpendicular wavenumber 
k„, the jump conditions are valid for each Fourier coefficient 
z^{x) - z^{x,kn) at x = and x - Lc, if and Lc are the two 
coordinates of the transition region. In the following the integer 
subscripts 1 and 3 will be used to label the layers as in Fig. 2, 
while the fourier modes will be labelled with the generic index 



2.2. Nonlinear model: Shell model for Reduced MHD 

In addition to the linear propagation of perturbations parallel to 
the loop mean field, we consider the waves to have perpendic- 
ular structure, so that the wave-vectors also have non-vanishing 
components in planes perpendicular to the mean magnetic field. 
In this transverse direction nonlinear interactions between differ- 
ent perpendicular modes occur, while the dynamics of the paral- 
lel propagation (for a given perpendicular mode) remains purely 
linear. This model, known as Reduced MHD or RMHD (Strauss 
1976), is believed to be well adapted to situations with a large 
uniform axial field Bo compared to perturbation amplitudes and 
strong anisotropy in the sense that the scales perpendicular to the 
field are shorter than the length of the coronal loop (Rappazzo 
et al. 2007): 



dt 



^Jp ox ^ p 



(14) 
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where we have taken identical kinematic viscosity and resistiv- 
ity, the density is uniform in the direction orthogonal to the field 
and the total pressure gradient guarantees the incompressibility 
of the fields via the Poisson equation 

Vi(p^) = -V • (zl ■ Vzt). (15) 

A second approximation consists in transforming the perpen- 
dicular nonlinear couphngs by replacing them, at each point of 
the X coordinate mesh along the mean field direction, by a dy- 
namical system defined in Fourier space, which allows reach- 
ing a very high Reynolds number compared to genuine Reduced 
MHD. liiis is known as the Shell model for RMHD or hybrid 
Shell model (Nigro et al. 2005; Buchlin & Velli 2007). The 
Reynolds number gain can be quantified as follows. Assume 
K is the perpendicular resolution (ratio from largest to small- 
est scales. Assume also the parallel resolution scales as K^^^. 
When passing from the RMHD to Shell RMHD the number of 
degrees of freedom changes from K^/^Log2{K) ^ K^^^ 

(see below). The CPU time required to describe the same large 
scale evolution is proportional to this number multiplied by K. 
Conversely, the resolution reachable goes as the CPU time T as 
T^/" in the RMHD case and as T^'^ in the Shell RMHD case, 
thus passing from a resolution /To to a resolution A"^'^^. The same 
is true for the Reynolds number (which goes as a power of the 
resolution K), hence typically passing from 10^ to 10^. 

Coronal heating driven by photospheric motions has been 
studied using both RMHD and RMHD Shell models in a one- 
layer atmosphere (corona) version, with uniform Alfven speed 
and closed (line-tied) boundaries, i.e. imposing the photospheric 
perpendicular velocity at loop foot points. Here we will use an 
RMHD-Shell model, but in the three-layer context, that is, in- 
cluding the linear jump laws defined previously at the transition 
region for each of the perpendicular wave number. 

The Shell model is characterized by the number -n 1 of 
perpendicular wave modes, each being characterized by a per- 
pendicular wave number, with amplitudes (the direction of 
the wave vector is not specified in the model), with the follow- 
ing discretization: 

k„ = Tko n = 0...N (16) 

Starting from the RMHD equations, one can write the follow- 
ing simplified equations (see Buchlin & Velli 2007 for the full 
equations with non homogeneous density): 

dtz: + Vad^zt - r; - vklz: 

d,z~ - Vad^Zn = T- - vklz;, (17) 

where Va is either y° (chromosphere) or (corona), y is the 
kinematic viscosity (equal to the magnetic diffusivity), and the 
r* are the nonlinear terms which are a sum of terms of the form: 

- AkmZpZq with m, p, q being close to n (see Biskamp 1994; 
GiuUani & Carbone 1998 for the fuU expression of T*). 

From the basic Eqs. 17, one can deduce the (exact) energy 
budget equation of a flux tube of length L, section nf^^ and den- 
sity p, assumed constant, as: 

AEIAt = Fit) - D{t) (18) 

E is the total energy, F the energy flux and D the energy dissi- 
pation rate defined as: 

1 

E = M— Ax (m^ + b^lp) 
Jo 



= M^J dx [iz^f + iz-f] (19) 

F = MVa^ [iz^f - (ztf + (zlf - (Zof] (20) 

D ^ M— \ dx Yvklizt^+Zn^) (21) 

Here M - nf^^Lp is the mass of the loop system, and b^/p, 
iz'^)^, iz~)^ are the sum respectively of the energies per unit mass 
in all the modes n = 0..N. When applying Eq. 18 to the corona, 
we take pc - e^po, Va = and the subscripts (integration inter- 
val) 0, L represent the left and right coronal boundaries respec- 
tively (not including the chromosphere when it is present). Note 
that the parameter /j_o stands for the largest scale available in the 
simulation, which is always in all runs l±o = 41 ±. In the follow- 
ing, we will use the notations E,D,F as defined in Eqs. 20-21 
but always normalized by the total mass M of the loop system, so 
obtaining average energies and dissipation rates per unit mass. 

Several remarks are in order. First, the nonlinear terms don't 
appear in the energy budget equation 18, because the total en- 
ergy is conserved by nonlinear coupling, as well in the Reduced 
MHD equations as in the presently used Shell model version 
of the equations. Second the energy accumulated or lost by the 
corona is not directly controlled. Indeed, the energy flux entering 
the corona (Eq. 20, see also the more explicit Eq. 25 below) is 
determined by the difference between the incoming and outcom- 
ing energies at the two transition regions; as will be made clear 
in the next subsection, the boundary conditions fix the incoming 
amplitudes, possibly in terms of the outcoming amplitudes, but 
not the energies. 

2.3. Boundary and jump conditions for three- and one-layer 
model 

As a rule, boundary conditions are defined by imposing the value 
of at X = (the rightward propagating wave amplitude) and 
the value of z,7 at the boundary x - L (the leftward wave ampU- 
tude). The jump conditions at the transition region 

2.3.1. Closed model (line-tied) 

The loop contains only the corona. The usual closed or line-tied 
model has for boundary conditions: 

z:(x = 0,0 = 2f/,^(f)-z;(x = 0,f) 

z„(x ^L,t)^ 2Uj;(t) - zl{x = L, f) (22) 

The previous equation results from Eqs. 12 with a = 1, s 
z;^(x = 0, 0, and = z^(x = L, f). The z\ (z^) signal in the 
corona is obtained by prescribing the velocity amplitude 
(U^) of each mode n at the boundary x = (x = L). One checks 
from Eq. 18 that when mo,„ = 0, then the energy flux (Eq. 21) 
injected in the domain becomes indeed zero. 

2.3.2. One-layer model 

In this first model including leakage, the chromosphere is ex- 
cluded from the domain, the domain boundaries coinciding with 
the T.R., as in the closed model. The boundary conditions now 
take the wave jump conditions (a < 1 in Eqs. 12) expUcitly into 

account: 

z:(x = 0, = (1 + a)Z^it) - az„ix = 0, /) 
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z-ix = L,t) = il+ a)Z-{i) - azlix = L, t) (23) 

The quantities and now denote the prescribed wave am- 
plitudes entering from the chromospheric side of the transition 
region. With reference to Fig. 2 we have Z^ = z\ and Z„ = z'^. 



and is the only term that can cause leakage in the closed model. 
Note however that, in the closed case, it is non zero only for the 
injected modes (which are at large scales, see next subsection), 
due to the presence of the z\ factor: this strongly Umits the leak- 
age in the closed case. 



2.3.3. Three-layer model 



2.4. Parameters and timescales 



In this second model allowing leakage, the chromosphere is re- 
ally included within the domain; in that case the boundary con- 
ditions (at the photosphere) are chosen to be purely open (a = 
in Eqs. 12 or equivalently in Eq. 23): 



ztix ■ 



■■ 0,t) 
L,t) 



■Kit) 
■Kit) 



(24) 



The boundaries are open in the sense that incoming waves are 
defined independently of outgoing waves, which in turn generate 
no incoming wave, so that they escape freely from the domain: 
perturbations coming from the loop reach the boundary and dis- 
appear below the boundary without reflection. Wave reflections 
and transmissions continuously occur within the domain at the 
location of the transition regions: there we apply the jump con- 
ditions (Eqs. 12), for each perpendicular mode n. 

The three-layer model and the one-layer model with partially 
reflecting boundaries are parametrized by the same number e, 
the photospheric/coronal Alfven speed ratio. The two models 
thus both include the transmission and reflection of waves by the 
transition region, but have an important difference. In the one- 
layer model, the chromospheric input is specified, as the T.R. 
coincides with the boundary of the domain. Instead, in the three- 
layer model, the chromospheric input {z\, Zj in Fig 2) is not 
prescribed, since the (prescribed) photospheric input has been 
modified by turbulence during its propagation through the chro- 
mosphere. Both models have specific advantages: the three-layer 
model has more internal degrees of freedom, as it shows two dis- 
tinct (but coupled) turbulent layers, one in the chromosphere, the 
other one in the corona; on the other hand, the one-layer model 
is more directly comparable to the closed line-tied model: the 
domain is the same (the corona), only the boundary conditions 
change. We will study both models, with some emphasis on the 
three-layer model. 

In the simulations we present, forcing is applied by injecting 
upward propagating waves only at the left loop foot point; more 
precisely, the leftward propagating amplitude at the right pho- 
tospheric foot point Z^(f) will be maintained zero in Eq. 23-24. 
In the particular case of the closed model (Eq. 22), this means 
that the velocity at the right foot point U^{t) was kept zero. In 
this case, it interesting to write down the expression for the net 
coronal energy flux : 



F oc SJ (1 + af\z\\- - 2a(l + a)Re(zt • zl*) 



-(l-a%£|2-(l-a%^ 



+ |2 



(25) 



In the previous formula, the ★ denotes the complex conjugate, 

indices n are assumed for each variable; we have used the no- 
tations of Fig. 2, in order for the formula to apply to the three 
models. We consider in turn the diflerent terms in the right-hand 
side. The last two terms are always negative: they thus represent 
a pure leakage (and they indeed vanish for a = 1 , in the closed or 
line-tied model). The first term is always positive and represents 
the continuous energy injection. The second term is fluctuating 



The parameters of the model are the length of the chromo- 
spheric and coronal part of the loop Lch, Lc respectively, the 
photospheric-chromospheric Alfven speed, V^, the Alfven speed 
contrast e, the width of the loop Z^o, the turbulent correlation 
scale and the amplitude of the forcing at the left photosphere, 
Uq. In all the models we will always force by injecting an Alfven 
wave: t/o is the wave ampUtude which, generally, is not directly 
related to the photospheric velocity shear. Only in the closed 
model the two quantity coincide (see section 2.1). The input 
photospheric spectrum will be distributed on the perpendicular 
scales l±/2, l±/4, and will have a correlation time given by 
Tf, which completes the set of parameters. 
For all the simulation we have set = 700 m/s, Lch = 2 Mm 
(so that L scales with L, only), i.e. we assume that photospheric 
values are independent of the loop length and that all loops have 
a transition region. We will also set l± = /j_o/4 and Tf = oo. The 
rest of the parameters l±, L^, Uq, e define the following phys- 
ical time scales (i.e., input of the model): the leakage time, the 
coronal Alfven time, and the input nonlinear time (which rules 
the strength of the turbulence resulting from the driving): 



ti 
f 



IJilnUo) 



£tL 



(26) 
(27) 
(28) 



In this work we will also fix = 6 Mm and e « 0.02, thus only 
fj^^ will be varied at fixed f^j and t^, by changing the parameters 
Uq and Zj. (in a following paper we will study the effects of vary- 
ing the leakage and the Alfven time scales). From these charac- 
teristic times we define the following dimensionless parameters 
which measure the nonlinear term vs the two main linear effects, 
the Alfven wave propagation and the leakage: 



XL = tjtli^ 

XO = tfalt%L = ^Xl 



(29) 
(30) 



The parameter has been used by Dmitruk et al. (2003); 
Rappazzo et al. (2008); Nigro et al. (2008) to quantify the turbu- 
lent behavior in their studies of turbulence forcing with closed 
boundaries (corresponding to = oo). 

Fig. 3 shows the plane with;(fo in abscissa and;^'/, in ordinate. 
This plane is divided in four quadrants by the lines Xd - 1 ^ind 
Xl = ^- There are actually only three subsets left, as only the 
subset with e < 1, visible as the non shaded region of Fig. 3, is 
permitted, due to the stratification. Turbulence will be said weak 
(in the left part) or strong (right part), depending on XQ being 
smaller or larger than unity. In the two upper quadrants, which 
occupy most of the domain, leakage should be negligible. Only 
in the small (left) bottom region should leakage dominate over 
turbulent loss. 

The two curves (solid and dashed) represent each a family of 
coronal loops of varying length L, build from a two-temperature 
hydrostatic model (see appendix A) which leads to a function 
y^(L). The solid (resp. dashed) curve corresponds to a coronal 
temperature of 3 (resp. 0.8) MK, and the loop length L increases 
from bottom to top (i.e., with increasing ;jf/,) from 3 to 700 Mm. 
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10"^ 10"^ 10"' 10° 10' 10^ 

X0-'^*NL 

Fig. 3. Characteristics of typical solar loops compared with simulation 
parameters: xl versus xa- The turnover time is fixed to f^^ = 1000 s and 
we assume a two-temperature hydrostatic loop model (see appendix). 
The black solid and dashed lines are for the two coronal temperatures 
of 3 MK and 0.8 MK respectively. The circles are for the 3-layer runs 
A-/ (see Table 1). The diagonal dotted lines are the curve = {\le)xa 
for three values of e. 



As shown by dotted lines, most of the hot loops show an Alfven 
speed contrast of e = 0.003, about ten times smaller than the e 
value chosen for the 3-layer simulations (see circles). The choice 
of relatively large e values for the simulations comes from the 
requirement of having a reasonable value for the ratio of integra- 
tion time to single time step. 

Let us now give a brief account of the physical and numer- 
ical time scales. Taking for instance /j_ = 2 Mm, with N = 20 
perpendicular wave modes, the largest available perpendicular 
wavenumber will be k,„a_^ = 1 .600 1 /km. The shortest nonlin- 
ear time (evaluated at the maximum perpendicular wavenumber 
in the corona) will be, if Zc is the typical coronal amplitude (z 
denoting either or z") 

TnL = l/ikmaxZc) - eliKiaxUa) (31) 

where we have taken the resonant linear case (see next section) 
for which the wave amplitudes is larger by a factor 1 /e in the 
corona. Replacing by previous values and assuming e = 0.01 we 
obtain for the smallest nonlinear time: 

T^z. - 6 10-^ s (32) 

As a matter of comparison, we use N\\ = 10"* grid points to de- 
scribe space along the loop, so that, for a typical loop length 
L - 6 Mm, we obtain for the smallest linear time for parallel 
propagation in the corona: 

T|| = l/ik'f'V'A) = eL/(7rNn\i) - 2 10-^ s (33) 

Hence, the constraint on the time step comes from the perpendic- 
ular nonlinear time. Finally, at least in the linear case (see next 
section), the characteristic time for large scale evolution is the 
long leakage time 

fi = L/(V^) - 9 lO^s (34) 

Comparing Eqs. 32-34, we see that « 10^ time steps of a dynam- 
ical system with 2 x 20 x 10"^ degrees of freedom are necessary 
to achieve one (anticipated) characteristic evolution time of the 
system. 



Table 1. Parameters for the the simulations 
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Notes. The three subpanels refers to the 3-layer, 1 -layer, and closed 
models (from top to bottom respectively) in which boundary conditions 
are open (3-layer), half reflecting (1-layer ), and line-tied (close). For 
all runs: Tf = oo, = 0.7 km/s, L^/, = 2 Mm (except for runs H, /, L 
that have L^, = 1 Mm), e denotes the ratio of photospheric over coronal 
Alfven speed (it plays no explicit role in the closed model). L,. is the 
length of the coronal part of the loop, /j^ is the perpendicular largest 
injection scale (the injected energy is distributed ar scales Z^, /x/2, and 
/l/4). Uq is the amplitude of the input wave.;^'o is the linear to nonlinear 
time ratio (Eq. 30),;^^l is the leakage to nonlinear time ratio (Eq. 29). 

3. Phenomenology 

3. 1 . Linear coronal trapping and leakage 

Let us first recall the linear result in the zero-frequency case, 
i.e. when forcing is time independent: a transverse perturba- 
tion (here, any perpendicular mode) is subjected to successive 
transmission-reflection at the two coronal bases, left and right. 
Since nonlinear interactions are ignored, all modes show the 
same evolution. As shown in Grappin et al. (2008) for a loop 
with smooth variation of the Alfven speed, the level of and z" 
grows progressively in the corona, in such a way as to achieve 
over a long time scale ti the asymptotic values 

2+ ^-z- ^bjyfp^ Uole (35) 

In other words, the asymptotic solutions are a uniform magnetic 
field amplitude everywhere along the loop at equipartition with 
the photospheric energy density, and, as well, a uniform veloc- 
ity fluctuation everywhere along the loop. The asymptotic state 
is thus the same as that would be achieved if the plasma were 
completely transparent to Alfven waves (e = 1): 

M = f/o, b ^ bo ^ Bi,Uo/V° (36) 

although this happens on the long time scale ?£ = L/V^ and 
not on the short Alfven coronal time f^. (We assimilate here and 
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in the following the coronal length to the total loop length L). 
Typically, if the shear ampUtude is Uq = 0.1 m/s, and the mean 
field Bo = 100 G, then the equilibrium magnetic field associated 
with the shear is the equipartition field, that is, bo - 14.5G. 

3.2. Resonant response 

Consider the simplest case where the frequency of the photo- 
spheric input is either zero or resonant (that is, equal to n/f^, 
with n an integer > 0). The coronal field perturbation induced by 
the photospheric field perturbation Uo = bo/ ^o grows linearly 
with time until it saturates at a finite value because of the two 
damping losses, the linear leakage (with time scale ti) and the 
nonlinear turbulent damping (with time scale to): 

d,b = BoUo/L - b/tL - b/tD = BoUo/L - b/t„ (37) 
= bo/tL-b/tr, (38) 

where bo - 14.5G is the photospheric magnetic perturbation, 
and tri is the effective damping time: 



tr, = (l/tL+l/tD) 



-1 



(39) 



Note that Eq. 38 we have rewritten the first term using the defi- 
nition ti^ = L/Vgin order to illustrate the fact that, in the absence 
of dissipation {to = oo, t,, = ti) the trapping and leakage times 
are equal. 

The stationary solution is for the coronal field perturbation: 



b = Uo {tr,Bo/L) 
= bo/il + tdtD) 



(40) 
(41) 



One sees that the coronal response is maximal (equal to the 
photospheric value bo - 14. 5G) when no turbulent damping is 
present (to » fe). In the other limit (ti » to), turbulent damping 
limits the coronal field to a fraction bo'.b ^ boto/tL = toBoUo/L. 

Relation 41 may be rephrased in terms of energy per unit 
mass as 



E = EoitJQ 



(42) 



with Eo = l/2(f/2 + bllAnpo) = f/^. In the case where to ^ 
tL, Eqs. 41-42 have been already given by HoUweg (1984); as 
remarked by Nigro et al. (2008), they are valid also for the zero 
frequency case (see also Grappin et al. 2008), the only difference 
being that in the latter case magnetic energy is dominant in the 
corona, while in the case of non zero resonance coronal magnetic 
and kinetic energies are at equipartition. 

A last remark concerns the use of Eqs. 38 and 41 (but 
not Eq. 42). Caution must be taken when applying the line- 
tied limit, ti - CO, since the trapping time, appearing as ti in 
these equations, is finite and fixed. Hence the explicit forms, 
Eqs. 37 and 40, are better suited to understand the difference 
between the opened and closed models. In particular one sees 
that the coronal magnetic field grows linearly with time in ab- 
sence of dissipation (Eq. 37) while, when dissipation is present, 
it can grow well beyond the leakage-limited value bo (Eq. 40), 
since the loss timescale t,, has no upper limit^ 

3.3. The general case 

In general, the signal injected into the corona is not necessarily 
resonanl and more generally not monochromatic. To quantify 



^ As we will see in the closed model the dissipation timescale can be 
larger than the nominal value 11 (see Fig. 13), thus leading to Z; > bo. 
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Fig. 4. Sketch of the linear coronal energy gain, log(£'/iio) as a function 
of frequency. £o is the input photospheric energy at each frequency. It 
is assumed that « t,,, and only the non-wkb portion of the spectrum 
is shown. The injected spectrum is also indicated as a shaded areas in 
arbitrary scale, for the case of flat, intermediate, and line spectrum (in 
increasing gray scale order). If leakage is neglected — > oo) dissipa- 
tion determines the height and width of resonances, since it is the only 
mechanism that limits the energy accumulation. For weak dissipation 
(and generally for — > oo) the resonances tend to delta functions. 



both the trapped energy and its dissipation rate we need to know 
how the time-dependent energy input is distributed between res- 
onant and non-resonant frequencies. We thus introduce the cor- 
relation time of the energy entering the corona or equiv- 
alently the width of the injection spectrum l/T^gr' which is a 
priori unknown. 

Recalling that the resonant lines are spaced each 1 /f^, each 
with a width equal to the inverse of the damping time tj^, then we 
may distinguish several cases depending on the portion of the 
excited spectrum (see Fig. 4): 

- Flat spectrum : T^^^ < f^. Since negligible energy is transmit- 
ted outside the resonant lines (enlarged by damping) com- 
pared to the energy transmitted for frequencies within the 

lines (anti-resonances), this leads to a filling factor equal to 
f^ltjf compared to a spectrum made of only resonant frequen- 
cies (Eq. 42). 



- Intermediate r°„, : < 



< t,j. Then the filling factor 



is T^or/^n as only the zero-frequency resonance and the first 
anti-resonance are excited. 
- Long correlation time or resonant spectrum: t,, < T^^^. This 
coincides with the linear resonant gain Eq. 36 if tL « to- 



Finally: 

E = Eo{t„/0 
E = EoitJO Tl,lf„ 
E = Eo(t„/f,f 



\^ cor ^ 'a/ 

(fa < TL < tn) 

it, < T'cor) 



(43) 
(44) 
(45) 



Note that we transformed Eq. 44 which originally reads E = 
Eoittflfaf' T^^„Jtri. Equations 43-44 have been derived for negli- 
gible leakage (f^, - to), in the strong turbulence case by Hollweg 
(1984) and in the weak turbulent case by Nigro et al. (2008). The 
relations proposed here extend these early findings by including 
the case where leakage dominates turbulence, and the case of 
very weak turbulence (resonant spectrum). 

To make these expressions explicit, one should express the 
unknown parameters in terms of control parameters. It is tempt- 
ing for instance to identify r',',^ with f*^^: then the three regimes 
correspond respectively to strong turbulence (xo > 1)) weak tur- 
bulence (xo < 1)> and weak dissipation (xl < 1)- 
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li. to 
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Fig. 5. From left to right: Run Dd (closed model), Run Dil (1-layer model) and Run D (3-layer model). For all the runs xo ~ 0.2, for the 1-layer 
and 3-layer model « 11.. Top panels: time evolution of the net energy flux F. Bottom panels: time evolution of the coronal energy (E) and 
dissipation below (D). Time is normaUzed to the input nonlinear timescale f^^, energy is normaUzed to the injection energy at the left boundary 
Etr = U^, dissipation and flux are normalized to HoUweg expression Dg = Ejitlf^. 



Indeed, the latter condition can only be satisfied if ti < to 
if we exclude the possibility that to < t^^. Thus assuming 
^cor - ^NL' '•^^ spectrum coincides with the linear reso- 
nant gain, and can only be reached by imposing xl < 1 • If 
the injection spectrum has a finite width, a possible choice is 
= min(r/,fO^) as suggested by Malara et al. 2010. If 
Tf > t^^ we fall in the previous case. If instead Tf < f^^, the or- 
dering considered in Nigro et al. (2008); Malara et al. (2010), the 
Une spectrum is not achievable. However, as we will see, the cor- 
relation time may also be given by other timescales, like the leak- 
age time ti and/or the chromospheric crossing time f^'' - LchlV^^. 
The remaining (difficult) task is to express the dissipation time 
to (since = min{tD, ti)), in terms ofxo ^^Xl via the coronal 
nonUnear time. We wiU come back to this point later on. 

3.4. Dissipation in ttie strong and weak regimes 

In the strong turbulence case (xo > 1) dissipation is expected 
to dominate on leakage and a simple exphcit expression of the 
dissipation rate is obtained after replacing = to in Eq. 43 
(HoUweg 1984): 

D = E/tD = £o/f« (46) 

This relation is attractive, as it leads to a universal result: the 
heating rate per unit mass does not depend on the detail of turbu- 
lent dissipation, since it depends only on the length of the loop 
and the photospheric energy. However, this universality is lost 
when we turn to the weak turbulent regime, Xo < Eq. 44, 
which we have seen is probably prevalent in the corona (Fig. 3). 
To extrapolate the previous expression (Eq. 46) to the weak 
regime with xo < 1. we identify T^g^ with in Eq. 44 and 
still adopt to < ti- 

D = EltD^Edtl{\lx^) (47) 

This predicts that, the weaker the turbulence regime, the higher 
the dissipation. We will see that both relations 46-47 are reason- 
ably satisfied if we use the line-tied limit, but not in the more 
realistic open case. In the open case, we will find that actually 



HoUweg's expression (Eq. 46) holds more or less both forxo > 1 
and;t'o < 1, which requires to admit that to > ti in the weak 
regime, i.e., that the dissipation time becomes very long as tur- 
bulence becomes weaker. 



4. Results 

In the following we will compare first the different models in a 
weak turbulence case, the most probable for coronal conditions. 
Then we will focus on the 3-layer model, comparing the weak 
and strong turbulence regimes. 

4.1. How leakage ctianges turbulence: the weak turbulence 
case 

We consider here a weak turbulent case with xo ~ 0.2, and 
compare the closed, 1-layer, and 3-layer models. The runs are 
Dci> DiL, andD respectively in Table 1; in the open models 
Xl ~ 1 1 so we expect that turbulence is the main channel for 
energy loss in all models. Because of this, we should not expect 
significant differences between the closed and the 1-layer run. 
However, we might perhaps find differences due to the different 
forcing (from now on we will use forcing to mean injection into 
the corona) between the 1-layer and the 3-layer runs, recalling 
that forcing is constant in the first case, and time-dependent in 
the second, due to the possibility of a chromospheric turbulence. 

The time evolution of the corona in the three models (from 
left to right) is summarized in Fig. 5 where the entering energy 
flux F (top panel), the total energy E and dissipation D (bottom 
panel) are shown (see Eqs. 19-21). Time is normalized to the 
input nonUnear timescale, energy is normaUzed to the input 
coronal energy Etr = ^^1^, the dissipation and the flux are nor- 
malized with respect to Hollweg expression, Dh = Etr/ fa (for 
the closed and 1-Iayer model ztr = Ut), for the 3-layer model 
Ztr is the measured quantity z\ that is not directly controlled by 
the boundary conditions). 

A quick look at the energy flux curves shows a sharp con- 
trast between the closed run and the open 1-layer run. While in 
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Fig. 6. Coronal kinetic energy spectrum (dashed line) and magnetic en- 
ergy spectrum (solid line) for runs D^, Dil, and D. Wavenumbers are 
in units of 1/Mm, in the top x-axis the corresponding shell numbers are 
indicated. The spectra are averaged in time and space, the normalization 
is in arbitrary units (spectra are also rescaled) . The symbols on the E/, 
spectra indicate the first shell number for which t^j^ik) = \jku{k) < f^. 



the closed case, the coronal energy flux is almost always pos- 
itive, in the open case, it is constantly oscillating around zero, 
although with a positive mean value flux. This has an imme- 
diate corollary: the energy level shows much lower values in 
the open case. Another corollary is that the dissipation rate it- 
self, i.e., coronal heating, is reduced by a factor ten. This ten- 
dency is sharply enhanced in the case of the three-layer model 
which shows a further reduction of a factor 5. Another remark- 
able difference appears in the 3-layer model, which accounts for 
the chromospheric turbulence. The energy and the energy flux 
display quasi-periodic oscillations that are absent in the closed 
and 1-layer model whose energy time-series are shaped by the 
time-independent forcing. Such oscillations have a periodicity of 
the order of one leakage time or smaller (see the top horizontal 
axis in the bottom panel). However we cannot rule out that their 
origin lies in the chromospheric turbulence. Indeed, the period- 
icity happens to be close to two chromospheric crossing times 
2fi^ = 2Lch/LctL - 2/3fL, which we interpret as the timescale 
necessary for waves injected from the left footpoint to leave the 
chromospheric layer (a round trip of the chromosphere). Most 
probably such oscillations are due to the coupling of the chro- 
mospheric and coronal turbulence and both timescales matters, 
as we will see in section 4.3. 

Fig. 6 shows the time and space-averaged kinetic and mag- 
netic spectra in the corona for the three models. One sees that 
all cases show well-developed power-law ranges, plus a mag- 
netic hump at large scales. The (common) forcing range is rep- 
resented by a gray vertical band. Note also the symbol on the 
magnetic spectrum which marks the large scale for which the ef- 
fective nonlinear time computed on the rms velocity at that scale 
is smaller than the Alfven crossing time. The only significant 
difference visible between the three spectra is that the magnetic 
peak is located at the largest forcing scale for the open runs, 
while it has migrated to a scale larger by a factor two for the 
closed run. This indicates that an inverse transfer is active in all 
cases, but that it is more active in the closed case, or also possi- 
bly that it has been hindered by leakage of the largest scales in 
the open cases. 



Fig. 7. Run D (3-layer open model, weak turbulence case): spatial distri- 
butions of fluctuations (top and mid panels) and turbulent heating (bot- 
tom panel). The time averaged rms amplitude (in km/s) are plotted as 
a function of loop coordinate (in Mm) for z*, z (top panel, solid and 
dashed line respectively) and for hj -^^npc, u (bottom panel, solid and 
dashed line respectively). The time-averaged heating rate (in arbitrary 
units) is plotted in the bottom panel as a solid Une, also a snapshot is 
shown in dotted Une. 



We are thus forced to conclude that, in the open models, de- 
spite the fact that;jfL > 1, the energy accumulation is limited by 
leakage. This means that the nonlinear timescale is a sharp 
under-estimation of the real dissipation timescale. We will come 
back on this point in the following. 

4.2. The 3-layer model: chromosphere vs corona 

We detail here the structure of the open 3-layer model in the 
weak turbulence case. In particular we compare the chromo- 
sphere and corona. We show in Fig. 7 the spatial profiles in the 
corona and chromosphere of the fluctuations and of the dissipa- 
tion rate. The top panel shows the time average of the rms value 
and z~ amplitudes with z^ms defined as: 



The mid panel shows the time averaged rms values of velocity 
u and magnetic field in km/s units (b/ y[p). The bottom panel 
shows the time average (solid line) and a snapshot (dotted line) 
of the heating rate. 

In spite of the presence of turbulence (as revealed by the 
spectra examined above), the rms amplitudes of all quantities 
are seen to be remarkably smooth functions of loop coordinates 
except of course at the T.R.. Main features are : (1) the mag- 
netic field amplitude in the corona and chromosphere are ac- 
tually comparable (the magnetic field ampUtude plotted in the 
figure is bj -^^^np, hence a factor of about 1/e = 50 between the 
coronal and chromospheric values); (2) the velocity contrast is 
significantly larger than unity but much smaller that the mag- 
netic contrast (in units of velocity), and its coronal profile has a 
simple form (3) the and z" levels are comparable in the corona 

(1) implies that the main part of the magnetic energy trapped 
in the corona is actually close to the linear state of zero resonance 
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Fig. 8. Run D. Contour plot of the spectra E-(x,k„) (snapshots) for z* 
and z" (top and bottom panels respectively) compensated by k^'^ in the 
chromosphere (left panels) and in the corona (right panels). Ordinate: 
shell number = log2(^„/fco)- Abscissa: coordinate x along the loop 
in Mm. The contours have different ranges in the chromospheric and 
coronal layers to better highlight their structures. 



(in the linear case with zero frequency the asymptotic coronal 
magnetic field fluctuations is equal to the photospheric field, see 
Grappin et al 2008). (2) the coronal profile of the velocity field 
is actually close to the profile of the first linear resonance (Nigro 
et al. 2008). (3) allows full nonlinear coupling which is compat- 
ible with the existence of a developed spectrum. 

Nigro et al. (2008) had already found in the closed case that 
the characteristic linear resonance profiles of the coronal cavity 
are not deeply affected by the presence of a nonlinear cascade. 
It appears that the same linear resonance profiles are'nt affected 
by leakage either. 

Finally, the time-averaged profile of the average dissipation 
rate per unit mass (bottom panel in Fig. 7) shows that the chro- 
mospheric dissipation remains negligible, and also that the left 
T.R. (i.e. above the foot point where energy is injected) is dissi- 
pating at a slightly higher rate than the other foot point. A typical 
snapshot (in dotted line) is also shown, providing a hint of the 
substantial intermittency of the heating rate, both in space and 
time. 

The turbulent activity of both the chromosphere and corona 
are shown in figure 8, in which we plot snapshots of the and 
spectra E'^ix, k±) (top and bottom panels) respectively in the left 
chromosphere and corona (left and right panels). The spectra are 
compensated by k^^^ in both layers. Note that the motivation for 
plotting spectra instead of u and b spectra is to make clear the 
respective contributions of the chromosphere and corona to the 
spectral formation, as the directions of propagation are identifi- 
able for 2=^, not for u and b. Note that only the left chromosphere 
has been represented (since in the right chromosphere the evolu- 
tion is purely linear, due to the absence of z" input from the right 
foot point), that its length has been enlarged to make its structure 
more conspicuous, and that the contours have different ranges in 
the chromopshere and in the corona. 

One can see in the figure something like the trajectory of tur- 
bulence from the left foot point to the corona and all the way 
back (so, one begins from the top left panel and proceeds clock- 
wise). First, in the chromosphere the onset of turbulence does 
not take place immediately starting from the left foot point : the 



spectrum (top left) first shows only the 3 injected scales (seen 
as a red-yellow band), and only very progressively adds smaller 
scales (first seen as a blue haze). Spatial intermittency then ap- 
pears about in the middle of the chromosphere in the form of 
small scale filamentary structures. Note that this coiTesponds to 
a travel time AL/V" ^ 1100 s, which is close to the nonlinear 
time f^^ ^ 800 s. 

In the corona (top panels) one sees on the contrary no large 
parallel gradients: as was seen previously with the rms and 
z^ energies. A conspicuous feature of the coronal spectrum is 
the hump appearing as a red ribbon which is displaced towards 
large scales (when compared to the peak in the chromospheric 
injected spectrum, top left). This again reveals the inverse trans- 
fer already noted above in Fig. 6. 

Finally, one sees in the bottom left panel that the wave leak- 
ing from the corona makes the z" chromospheric spectrum look 
much more developed than its z^ counterpart. 

4.3. The 3-layer model: Increasing turbulence 

We now increase in the 3-layer model the turbulence strength 
Xo from 0.04 to 4.8 (runs A, D, F, H). This is achieved by 
decreasing the nonlinear time, while the leakage time is fixed. 
Hence even though we have already seen that the nonlinear time 
is clearly a strong lower bound for the dissipative time, again, 
one should expect the open model to match at some point the 
closed model in the limit » 1. This point will be considered 
again in the discussion where the properties of all models are 
summarized. 

In Fig. 9 we illustrate how the dynamics changes in the open 
3-layer model when increasing ;t'o- The left panel shows the rms 
magnetic field amplitude in the corona normalized to b^, the lin- 
ear zero frequency solution, while the two other panels show the 
(space and time averaged) total energy spectra respectively in 
the chromosphere and the corona. 

The main points are (1) When the nonlinear time is too 
large (very small xQi run A), one sees that turbulence has no 
time to develop before reaching the corona. Both the chro- 
mospheric and coronal spectra remain largely devoid of small 
scales. Dissipation is thus negligible. The asymptotic level of 
the magnetic field is close to its 15 G linear value, the growth of 
brms being extremely regular and devoid of any small scale fluc- 
tuations. All this happens in a leakage time. (2) Decreasing the 
nonlinear time progressively decreases the asymptotic coronal 
field. Its growth becomes now chaotic, the signal in the left panel 
showing a whole spectrum of frequencies, with most remarkably 
periods close to the leakage time for the two intermediate values 
of Xq^ but also periods close to two chromospheric crossing time 
IfJ' - l/3tL for the strongest ;^fo (run H in the left panel, see 
for example the range t/ti G [3,4]). (3) At reasonably large ;t'o, 
the coronal spectra are developed. However, the chromospheric 
spectra are significantly steeper. In the chromosphere, the slope 
is close to 1.8, while it is close to 1.7 in the corona. (4) Note also 
that the chromospheric spectra are devoid of the humps which 
appear in the coronal spectra. 

We thus conclude that too weak a cascade does not change 
linear zero frequency results at all, and that there is axo threshold 
above which turbulence has common properties. There are slight 
differences in the chromosphere and corona, the main ones being 
the large scale coronal peak, and a slightly different slope. 

We now examine frequency spectra. We have computed fre- 
quency spectra of izf,„s)' at each position along the loop; then 
they have been separately averaged in the corona and in the chro- 
mosphere. The original time series has been windowed with the 
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Fig. 9. Runs A, D, F, and H with increasing ^o- Left panel: growth of rms coronal magnetic field (normalized to its asymptotic linear value 
bg = 15 G). Middle and right panels: total energy spectrum E{k) averaged in time and space in the chromosphere and in the corona respectively. 
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Fig. 10. Runs D and H (weak and strong ;i'o): Frequency energy spec- 
tra E^if), E^(f ), computed by taking the Fourier transform of the 
(z*„j(f))" at each plane and then space averaging separately in the chro- 
mosphere (bottom black and red lines) and corona (top black line). The 
E~ (not plotted) and E* coronal spectra are indistinguishable. The verti- 
cal lines mark some relevant timescales: dotted lines for the resonances 
at n > 1, dotted-dashed lines for the (round-trip) chromospheric cross- 
ing time, 2fJ', and long-dashed lines for the leakage time, t^. 



harming procedure and the zero frequency is also displayed as 
the lowest frequency in the plot (Fig. 10). The coronal and 
z" have practically the same spectrum, so we plot only in the 
corona, as well as the chromospheric spectra of and z" ■ 

One sees on the coronal spectra the appearance of peaks 
close to (but not coinciding exactly with) the resonant frequen- 
cies: n/tg, these harmonics being marked as dotted lines. This 
confirms again that the quasi-linear trapping properties are not 
strongly affected by leakage. 

For weak turbulence (left panel) the spectra are dominated 
by the lowest frequencies, the input zero-frequency and a low- 
frequency bump. As the strength of turbulence is increased (right 
panel), more energy goes into finite-frequency resonances, some 



of them becoming as energetic as the low frequency part of the 
spectrum. Note that the location of the low-frequency bumb cor- 
responds roughly to two characteristic timescales, the leakage 
time, ti, and two chromospheric crossing time, It^^'', which we 
interpret as the signature of the turbulence activity in the coro- 
nal and chromospheric layer respectively. In run D (low turbu- 
lence) two distinct bumps appear in the chromopheric spectrum 
E^i^ at frequencies l/ti and 1/2 fj\ while in run H (strong tur- 
bulence) the bump lies in between them. In the coronal spectra 
(and also in E'^j^) the bump is somewhat wider, possibly show- 
ing a coupling with the (zero and finite frequency) resonances. 
The importance of both timescales points out the fact that the 
low-frequency spectrum in the 3-layer model is affected by the 
coupling between the coronal and chromospheric turbulence. 



5. Discussion 

We have studied in this work the problem of heating coronal 
loops by forcing a photospheric shear through the injection of 
Alfven waves at one foot point, and waiting for the injected 
Alfven waves to be transmitted into the corona, to be trapped 
(amplified), and finally to dissipate due to turbulence. We used 
for this a simplification of the RMHD equations that exploits 
Shell models to account for the perpendicular nonlinear cou- 
pling. 

In contrast to previous work, we have considered a finite 
Alfven speed difference between the photosphere and the corona 
and freely propagating waves deep in the photosphere, thus al- 
lowing energy to leak back to the chromosphere and to deeper 
layers of the solar atmosphere. We found that, although leak- 
age doesn't change dramatically either the quasi-linear trapping 
properties of the corona or the spectral properties of turbulence, 
it alters strongly the level of the energy trapped in the corona and 
the resulting dissipation. 

We have seen that the coronal energy and dissipation rate 
vary, depending (i) on the ratio xo of the linear Alfven cross- 
ing time by the nonlinear input time fj^^ (ii) on the choice of the 
atmosphere model, allowing or not the development of a chro- 
mospheric turbulence. We now present systematically these vari- 
ations and discuss them. We will present on an equal footing 
the closed model, 1 -layer and 3-layer models. The closed model 
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Fig. 11. Coronal energy per unit mass normalized by the T.R. injection 
energy as a function of the time ratio xo- Left: closed model results 
(symbols); Mid panel: 1-layer model results (symbols); Right: 3-layer 
model results (symbols); dashed line: resonant scaling; dotted line: in- 
termediate correlation time scaling; solid line: small correlation scaling. 
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Fig. 12. Dissipation normalized by Eq. 46 vs xo for the three models. 
The vertical dotted lines mark the values = e (f^ = tt] and;^;'o = 1. 
The dashed line is the prediction from Eq. 5 1 for the 3-layer model. 



results should allow direct comparison with earlier work, and 
the comparison between the 1-layer and 3-layer models should 
make clear the effect of chromospheric turbulence. However, we 
consider that the 3-layer model is the more realistic of the three 
models. 

Note also that we will vary;t'() while fixing the Alfven speed 
contrast to e = 0.02 and the loop legnth to = 6 Mm. The 
dependence on e and Lc will be postponed to a further study. 

5.1. Energy 

We show in Fig. 1 1 how the energy trapped in the corona de- 
pends on the turbulence strength ;^'o in the different models (left: 
closed model, mid panel: 1-layer, right: 3-layer). In each case, 
we have shown how the results are fitted by the generalization of 
the Hollweg-NMV model as given by the three different possible 
regimes (Eqs. 43-45). For this purpose we normalize the energy 
to the coronal input energy which coincides with the input 
energy Ul in the closed and 1-layer model while it is a measured 
quantity for the 3-layer model (Zj^, = in Fig. 2). 

The regime is identified by the choice of the correlation time 
of the input spectrum, which is different in the different models. 
For instance, a constant signal in the closed model is a particu- 
lar case of the resonant regime (dashed line, line spectrum). The 
flat spectrum should be found in the case of a very short correla- 
tion time (solid line), and the intermediate regime by the dotted 
line. It is striking that the best fit, although largely imperfect, 
is always obtained by the resonant expression (zero-frequency 
spectrum), whatever the turbulence strength (i.e. both for small 
and large ;t'o)- 

The behavior of the system is therefore dominated by the 
n = resonance. This could be somewhat expected for the 
closed and 1-layer models, since the input signal is time inde- 
pendent, but less so for the 3-layer model, in which the chromo- 
spheric turbulence modifies the input frequency spectrum to the 
corona, at least for large ;^'(). Hence, it seems that the dominance 
of the low-frequencies is not caused by the particular forcing 
chosen here, but it is a consequence of the coronal activity itself. 



As we have seen, the coronal spectrum shows a bump at large 
perpendicular scale, containing very low frequencies; this bump 
is not present in the chromospheric spectrum (see fig 9). 

A last point concerns the model differences. While all three 
models attain the same energy level for xo > 1, it is seen 
that they strongly differ in the weak turbulence regime. In the 
closed model the energy grows until turbulence becomes effi- 
cient enough to balance the input energy. In the opened models, 
instead, the energy accumulation is prevented by the leakage and 
the resulting level is much lower. 

5.2. Heating rate 

The heating rate is shown in a single plot for all three models 
in Fig. 12. We show the dissipation per unit mass normalized to 
Hollweg expression, D/, - Ejjilf^, (Eq. 46) vs;^'(): diamonds are 
for the closed model, stars for the 1-layer model, and circles for 
the 3-layer model. 

In the weak turbulent part of the diagram ixo < 1), we find 
again the same ordering of the models observed above for the 
energy: a line-tied model leads to a heating rate inversely pro- 
portional to the turbulence strength (our results follow well the 
fit D oc {u^IQxq^'^ proposed by Dmitruk & Gomez 1999) while 
for the open models the heating rate goes down proportionally 
to the turbulence strength. We note however the sudden drop 
of the dissipation rate at very low xo (roughly corresponding to 
^ their equality being indicated by the dotted vertical 
line at xa - £), due to the absence of formation of a high- wave 
number spectrum in the open models: the coronal energy level 
is too low to trigger a cascade before fluctuations leak out of the 
corona. 

In the strong turbulent regime (xo > 1) the dissipation is 
about independent of the turbulence strength: the level of this 
plateau is common to the closed and 1-layer models which differ 
only in their boundary conditions, but it is lower for the 3-layer 
model. This is to be attributed to the chromospheric turbulence 
that reduces the coronal input and hence the actual strength of 
the turbulence, for given ;t^o- 
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Fig. 13. Dissipation time vs xo for the three models; note the different 
normalization in the left and right axis. The vertical dotted lines mark 
the values ^0 = f UZ, = tt) andxo = 1- 



5.3. Dissipation time 

Finally we consider in Fig. 13 the dissipation time, to - E/D, 
normalized to the coronal crossing time t'^j. For the opened mod- 
els the right vertical axis also shows the dissipation time nor- 
malized to the leakage time, which makes sense since the ratio 
e between the Alfven and leakage time remains fixed in the data 
shown here. 

The few points of the open models with very low turbulence 
strength ;^'o show very high values of the dissipation time, due to 
the undeveloped turbulence, as already mentioned above. 

Leaving apart these points, one sees that increasing ;t'o (i-e-, 
decreasing the input nonlinear time), the dissipation time de- 
creases as expected, but that, most remarkably, it stops decreas- 
ing when xo reaches about unity. The plateau is comparable for 
the 1 -layer and closed models, corresponding to ^ lOf^ in 
the 1 -layer model, while it corresponds to ^ to in the 3-layer 
model. 

It would certainly be a progress, both from the theoretical 
and the practical viewpoint, to understand how the heating time 
and energy are related. The difficulty is that the usual relations 
here are modified by the existence of the large-scale hump in 
the spectrum: one cannot consider that there is a straightforward 
(direct) cascade from the large to the small scales, as the energy 
is clearly "blocked" at large scales. 

That this is so, can be easily verified by comparing the 
instantaneous dissipation time with the nonlinear time t^i - 
1 / (ku) (with k - kj_ being the largest forcing scale, and u the 
rms velocity) which usually rules the direct Kolmogorov cas- 
cade. However, a strikingly good result (shown in Fig. 14 as a 
red line, while the standard nonlinear time is shown in blue) is 
obtained when using a factor V^Ju to increase the nonlinear time 
l/iku): 



to - l/(ku)(V^Ju) 



(49) 



The fit is good enough and works for all values listed in 
Table 1 for the turbulence strength xo for the 3-layer model. 
However, the fit is more hazardous for the 1 -layer model, and 
simply does not work at all for the closed model. 

This simple law for the three-layer model deserves some 
comments. Note that the expression in Eq. 49 reminds one of 



10.0 




Fig. 14. Run D: Comparing the instantaneous dissipation time with the 
slow Iroshnikov-Kraichnan time scale (red) and the standard nonlin- 
ear time (blue, below). Note that the agreement between the measured 
dissipation time and IK time is as good as here for all values of the 
turbulence parameter ;^fo in the case of the three-layer model. 



the delayed cascade time predicted by the Iroshnikov-Kraichnan 
phenomenology. This phenomenology was meant to describe the 
(delayed) cascade of interacting Alfven waves with wavevectors 
not perpendicular to the mean field. The delaying effect was sup- 
posed to work all along the cascade, leading to a specific spec- 
trum, different from the Kolmogorov one (3/2 instead of 5/3). 
Here, the situation is different. Indeed, the coronal spectrum 
adopts a slope close to 5/3, not 3/2, and, most probably (but this 
study is postponed to a later work), the characteristic time scale 
which rules the cascade here is the strong turbulence time scale 
1 /{ku). However, the energy is dominated by the spectral hump 
at large scale which is not ruled by the fast time scale, but by 
the delayed time scale Eq. 49. This comes from the effect of the 
resonant trapped linear modes which act to deplete the nonlinear 
coupling terms, in the same way as in the Iroshnikov-Kraichnan 
phenomenology, although only at large scales. 

5.4. Conclusion 

We focus here on the 3-layer model which is by far the most real- 
istic model of the three we have studied. We have increased our 
knowledge concerning the physical mechanisms at work since 
we know now that a) the coupling between the chromospheric 
input and that of the coronal cavity is close to that of zero- 
frequency resonance b) leakage always plays a substantial role 
c) the time scale of dissipation is long, being the Kolmogorov 
time 1 /{ku) reduced by the u/V^ factor. 

We also have found that a rough prediction for the dissipa- 
tion rate is given by the classical expression: D ^ Dh - ^^/l^a 
(Eq. 46, see Fig. 12) where ztr is the amplitude of the input fluc- 
tuation at the T.R. level. This result is a bit paradoxical, since 
this relation has been first obtained by Hollweg, on the basis of 
assumptions which are not verified in our simulations: (i) a short 
correlation time for the chromospheric input (ii) negligible leak- 
age. In our simulations the conditions are completely different: 
(i) long correlation time (ii) substantial leakage. The solution of 
the paradox lies in the rough compensation of several effects: (i) 
the energy level is decreased due to leakage and to the chromo- 



13 



A. Verdini R. Grappin M. Velli: Turbulence and leakage in coronal heating 



spheric turbulence, but largely increased due to resonance (ii) 
the dissipation time is increased due to the large-scale hump. 

Although the above heating rate contains the unknown T.R. 
level of input fluctuations, it can be used as a predictive law if we 
identify the T.R. input value ztr with the (imposed) photospheric 
value Uq. This gives the simple classical result 



D = Dho = Ulltl 



(50) 



If we consider again the 3-layer results and plot the ratio D/Dho 
instead of the ratio D/Dh as in Fig. 12, it is interesting to note 
that the result is not basically changed, but nevertheless the de- 
viation from the horizontal (here D/Dhq = 1) is a bit reduced, 
being hmited to at most a factor 5. 

Can we use our new knowledge to improve the prediction 
of the dissipation rate beyond the approximate law D ^ Dh ? 
Unfortunately the answer is no, because of our poor knowledge 
concerning the relation between the known photospheric input 
Uo and the (largely unknown) chromospheric input ztr as well 
as the coronal velocity fluctuation level Uc. If we bypass this step 
by replacing the unknown quantities (ztr,Uc) by f/o, using the 
zero-frequency resonant expressions with dominant leakage for 
the coronal fluctuations (be = Uo/e, Uc = Uq, Eq. 45) and the 
dissipation time (Eq. 49), we obtain for the dissipation rate per 
unit mass: 



D = blltD - Ulle\\lt%^)(UdVl) 
= Dh (L/L) {UolVlf 



(51) 



When using this expression with the parameter values of the 3- 
layer model as given in Table 1, one obtains the dashed line in 
fig 12. This is clearly not an improvement of the simple relation 
D/Dh - 1 . It is actually much worse, by direct comparison with 
the numerical simulation results (the circles in fig 12) but also 
from a more general point of view, since dissipation is largely 
beheved to grow with Bq, while Eq. 51 predicts the reverse (D oc 

To progress, we must understand how to relate the chro- 
mospheric and coronal velocity level to the photospheric one. 
We should also explore how the heating rate depends on all pa- 
rameters, in particular the Alfven speed contrast e and the loop 
length Lc. Finally, we should investigate if the properties of pho- 
tospheric turbulence, in particular the correlation time, modifies 
or not the coronal reaction. 
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Appendix A: Atmospheric model 

In order to obtain the parameter for "realistic" coronal loop, 
we have to determine the Alfven crossing time, or in other words 
the relation e(L) for a given V^. We model the loop as a semicir- 
cular cy finder of radius R, subject to a constant vertical acceler- 
ation g = GMq/Rq. The loop has constant cross-section and is 
threated by a uniform magnetic field. For simplicity the loop is 
assumed to be isothermal in the chromosphere and in the corona, 
the two temperatures are related by the jump at the transition re- 
gion (T.R) 



ns) = To + i (Jc - To) 



tanh 



(A.1) 



where s g {Q,nR\ is the coordinate along the loop, 5^"* and 
6tr are the position and width of the two transition regions, 
sf^ = Rsirr\htr/R), sf, = L - sf,. We set flie T.R. height 
h,r - 2 Mm and its width to 6tr = 0.2 Mm. We will consider 
two coronal temperatures, Tc = 0.8 MK, 3 MK in order to con- 
sider short loop (reaching the low corona) and longer loops. We 
finally assign to the "base" parameters, magnetic field, number 
density, and temperature, the following values: Bq = 100 G, 
no = 10'^ cm-3 and Tq = 4500 K). 

The density profile along the loop is obtained by solving the 
equation for the static equiUbrium 



Idp IdT g 

= ; COS ns/L 

pas T as T 



(A.2) 



and the cos function accounts for the projection of gravity along 
the loop. Varying the loop length we will find short loops that 
don't reach the T.R. heights (^ < 2 Mm) and long loops, that 
indeed reach the corona. For the former the density is found by 
direct integration of the above equation, while for long loops the 
equations will be solved numerically. 

The static loop model, according to its temperature profile, de- 
fines a relation e(L) that can be estimated by considering the 
temperature jump at the T.R. as a discontinuity and calculating 
the density at the loop apex s = L/2. By integrating from the 
photosphere to the T.R. and from the T.R. to the corona one finds 



ln(p(r-/po) = -gh,r/To 

MPclptr+) = -giR-h,r)/Tc 



(A.3) 
(A.4) 



where p,r+ are the densities just below and above the transition 
region. Assuming that the T.R. is in pressure equilibrium, the 
density jump is given by ln(pir+/ptr-) = ^(To/Tc) so finally one 
gets: 



21ne = - 



R-h„ 

He 



Hn 



+ ln 



To 

Te 



Ho 



In 



To 



(A.5) 



where we have introduced the density scale heights in the chro- 
mosphere and corona. Ho = To/g ~ 0.27 Mm and He = Tdg ~ 
60 Mm respectively, and make use of the definition e - (pelpo?'- 
From Eq. A.5 one can see that for long loops, R < He (used in the 
last equality), the density contrast is determined almost entirely 
by the T.R. jump and is independent of the loop length, except 
for very long loops that span a density scale height in the corona 
(L » 200 Mm). By numerical integration of Eqs. A.1 -A.2 we 
obtain for a long loop e - 0.004, which is a factor two larger 
than the estimate based on Eq. A.5, the discrepancy arising from 
the fact that the T.R. is not in pressure equilibrium. 

To obtain the solid and dashed black lines in Fig. 3 we 
use the relation e{L) = V^/V^(L) as found from the numeri- 
cal integration for the two coronal temperatures Tc = 0.8 MK 
and = 3 MK (here V'JL) = mcLxAVa(L,x)]). The maxi- 
mal Alfven speed is shown in fig A.1 as a function of the loop 
length: increases monotonically and then experiences a sud- 
den jump at around L = 4.2 Mm. After that jump it decreases 
slightly and then increases again monotonically. In the first part 
(L < 4.2 Mm) loop are short enough to remain in the first isother- 
mal layer (the chromosphere), where the density scale height is 
small. The jump at L ^ 4.2 Mm is determined by the fact that 
the loop reaches the height of the T.R.. In this thin layer the den- 
sity scale height is very small, and density drop very quickly. 
Loops with length between « 4.2 Mm and 4.4 Mm don't pen- 
etrate into the corona, remaining in the T.R., thus the Alfven 
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Fig. A.l. Coronal Alfven speed as a function of the loop length L 
for the hydrostatic-two temperature model of coronal loop. Solid and 
dashed line correspond to coronal temperature of 3 MK and 0.8 MK 
respectively 



speed increases even more, reaching a local maximum. The next 
part of the profile is characteristic of loops that reach the second 
isothermal layer (the corona), where the density scale height is 
large. 
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